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1.  Introduction:  Statement  of  the  Problem 

In  the  computer  simulation  of  combat  over  a  specific  geographical  area,  the  terrain 
is  typically  represented  using  finite  sets  of  isolated  data  points  that  we  will  call  digitized 
terrain  altitudes  (DTAs).  These  DTAs  are  obtained  by  measurement,  either  direct  or  indi¬ 
rect,  and  are  generally  arranged  m  a  rectangular  grid.  A  100- meter  grid  spacing  has  been 
typical  for  many  combat  models,  although  the  current  trend  is  toward  somewhat  smaller 
grid  sizes,  especially  for  applications  in  the  area  of  distributed  interactive  simulation 
(DIS).  In  any  case,  when  a  terrain  altitude  is  needed  at  an  (x,  y)-coordinate  ii:  between 
grid  points,  it  must  be  approximated  using  the  surrounding  DTAs.  Moreover,  since  the 
DTAs  are  isolated  points,  there  is  no  direct  method  for  obtaining  terrain  slope  information 
even  at  the  DTA  coordinates. 

The  "inverse"  problem  is  to  construct  in  closed  form  a  "smooth"  surface 
z  =  T{x,  y)  (with  the  properties  that  T{x,  y)  is  continuous  everywhere  and  has  continuous 
partial  derivatives  y)  and  Ty(x,  y)  almost  everywhere  [i.e.,  everywhere  except  a  set 
of  measure  zero])  that  "matches"  the  DTAs  to  within  a  prescribed  tolerance,  S  (i.e„  for 
each  DTA  point  (x,  y,  d),  \T{x,  y)-d\<  S).  If  such  a  function  could  be  defined,  then  at 
any  point  on  the  battlefield  both  the  terrain  altitude  and  the  slope  of  the  terrain  in  any 
direction  could  be  obtained  simply  by  evaluating  T(x,y),  Tj,(x,y),  and  Ty(x,y)  at  the 
desired  (x,  y)-coordinate. 

The  benefits  of  such  a  terrain  description  in  combat  simulation  are  obvious.  All 
calculations  related  to  terrain,  from  line-of-sight  calculation  to  troop  movement,  could  be 
done  more  simply  and  at  any  desired  resolution.  In  this  report,  we  show  how  to  con¬ 
struct  a  function  r(x,  y)  with  all  of  the  above  properties. 

Once  T(x,y)  has  been  constructed,  it  is  possible  to  add  realistic  "micro-terrain"  lo 
the  surface.  We  indicate  how  this  can  be  accomplished.  The  critical  tool  in  this  process 
is  the  Variable  Resolution  Terrain  (VRT)  model. 

2.  The  VRT  Model 

The  VRT  model,  developed  at  the  U.S.  Army  Ballistic  Research  Laboratory  (now 
the  U.S.  Army  Research  Laboratory)  in  the  early  1990s,  is  a  model  that  represents  basic 
topography  as  a  cemtinuous  surface,  capable  of  being  viewed  at  any  desirexl  resolution.' 
This  surface  is  the  superposition  of  individual  terrain  features  or  "hills,"  each  of  which  is 

'  Wald.  Joseph  K.  and  Patterson,  Carolyn  J.  "A  Variable  Resolution  Terrain  Model  for 
Combat  Simulation."  BRL-TR-3374,  U.S.  Anny  Ballistic  Resciirch  Laboratory 
Abeidcai  Proving  Ground.  MD,  July  1992. 
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described  by  a  closed  form  mathematical  function  that  is  continuous  everywhere  and  has 
continuous  partial  derivatives  almost  everywhere.  While  examples  of  many  different 
types  of  terrain  have  been  ccmstructed  using  tlie  VRT  model,  they  are  all  "generic"  in  the 
sense  that  they  were  developed  from  basic  principles  and  do  not  represent  any  particular 
piece  of  Earth  topography.  For  generic  terrain  creatic«,  the  distribution  of  the  sizes  of  the 
hills  is  based  on  the  idea  of  self-similarity  (i.e.,  invariance  with  respect  to  scale)  which  is 
embodied  in  the  power  law 

D=^Ks~^,  (2.1) 

in  which  j  is  a  dimensionless  scale  factor  associated  with  a  hill,  AT  is  a  constant  that 
depends  on  terrain  type,  and  D  is  the  areal  density  of  hills  as  a  function  of  s.  The  choice 
of  exponent  in  this  power  law  ensures  self-similarity  in  the  density  of  terrain  features. 
Integration  of  this  power  law  produces  the  cumulative  distribution  function: 

/  =  5 

I  Kt-^dt  =  K(\l  -  Us).  (2.2) 

^  *  ^min 

In  practice,  there  is  a  range  of  scales,  for  which  the  power  law  holds.  Thus, 

to  build  each  hill,  we  need  only  draw  a  uniform  random  number  u  from  the  interval 
(0. 1  -  [-Smin !  -Jni'x]).  With  the  hill  scale  factor  being  r  =  5;nin  /  (1  “  «)•  The  location  of  the 
hill  is  randomly  chosen  in  the  desired  area.  The  complete  terrain  surface  is  defined  by  the 
superposition  of  all  of  the  hills,  i.e.. 


T(x,y)=  E/*(x,y). 
The  form  of  a  single  hill  function,  /*(x,  y),  is  given  by: 


AU,  y)  =  st  h^expi- 


where 


Qkix,  y)  =  a,(a:  - 1*)^  -  02^  -^*)(y  -  7t)  +  aj{y  -  ;7*)^ 


(2.3) 


(2.4a) 


(2.4h) 


^2  =  I  sin  22*, 


(2.4c) 

(2.4d) 
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and 


^3 


(2.4e) 


This  formulation  is  a  revision  of  the  original  hiL  function  definition.^  Vaiying  the  param¬ 
eters  I"*,  /i*,  p*,  £*,  A*,  <7*,  and  Si  produces  hills  in  a  variety  of  sizes  and  shapes. 

Figure  1  sho^vs  a  simple  "3D  wire  mesh"  of  a  VRT  surface  generated  by  this  procedure. 


3.  Constructing  T(x,y) 

The  first  step  in  the  process  of  "fitting"  a  smooth  surface  to  a  set  of  DTAs  consists 
of  detrending  the  DTAs  by  fitting  a  plane  through  them  and  subtracting  the  height  of  this 
plane  from  each  of  the  DTAs.  One  way  to  define  this  plane  (but  not  the  only  way)  is  via  a 
least-squares  procedure.  This  produces  "residual"  DTAs. 

Next  the  highest  residual  DTA  (jc,  y,  d)  is  found  (or  the  lowest,  if  the  largest  DTA 
happens  to  have  a  negative  value)  and  we  define  a  hill  whose  center  is  at  Oi*  near  (x,  y) 
and  whose  height  is  equal  to  or  close  to  d.  An  iterative  loop  is  used  to  find  the  best  val¬ 
ues  of  the  parameters  in  the  construction  of  this  hUl.  [Although,  we  use  equation  2.4  to 
define  the  hill,  any  function  with  the  flexibility  of  equation  2.4  arid  the  same  smoothness 
properties  could  be  used.]  Here,  the  best  hill  is  the  one  that  leaves  the  smallest  modified 
residual  DTAs  in  the  neighborhood  of  (x,  y)  when  the  hill  is  subtracted  from  the  residual 
DTAs.  This  neighborhood  is  user-definable  and  will,  in  general,  depend  upon  terrain 
type.  This  hill  is  subtracted  from  each  of  the  residual  DTAs  to  create  a  new  set  of  resid¬ 
ual  DTAs. 

We  now  proceed  to  find  the  new  highest  residual  DTA  and  repeat  the  above  hill  fit¬ 
ting  process  until  all  of  the  residual  DTAs  are  smaller  in  absolute  value  than  S.  The 
smooth  surface  Tix,  y)  is  then  defined  to  be  the  sum  of  the  fitting  plane  and  all  of  the  fit¬ 
ting  hills.  By  the  method  of  construction,  Tix,  y)  is  guaranteed  to  satisfy  all  of  the  condi¬ 
tions  of  the  inverse  problem. 

To  illustrate  the  terrain  fitting  process,  we  selected  a  1  kilometer  square  of  terrain 
located  near  Denver,  Colorado.  The  terrain  data  were  provided  by  the  U.S.  .^rmy  Topo¬ 
graphic  Engineering  CenU  r.-  The  grid  spacing  for  the  DTAs  is  5  meters,  resulting  in  a 
total  of  40,401  DTAs  in  the  set.  Figure  2  contains  a  contour  map  of  the  sample  DTAs. 

^op.  cii. 

’Ray  Ncrvelle,  U.S.  Anny  Topographic  Engineering  Center.  Private  communication, 

March,  1994, 
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Figure  1.  Sample  VRT  Surface. 
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We  proceeded  to  construct  three  smooth  surfaces,  which  we  denote  by  T^{x,  y),  y), 

and  ri(x,  y)  corresponding,  respectiv'ely,  to  ^  5,  <5  ~  2,  and  5  =  1.  lli?  contour  maps  of 

these  surfaces  are  shown  in  Figures  3,  4,  and  5. 


The  increasing  quality  of  the  fit  as  S  decreases  is  evident  from  the  contour  maps. 
However,  the  best  value  of  delta  for  the  fitting  process  need  not  always  be  the  smallest. 
Rather,  it  should  be  the  largest  value  consistent  with  the  user’s  intended  application. 
This  is  due  to  the  fact  that  as  S  decreases,  the  computational  burden,  both  in  constructing 
T (x,  y)  and  in  using  it,  increases.  Table  1  contains  the  computational  .  atistics  for  the 
construction  of  T^{x,  y),  TjCx,  y),  and  Tx{x,y).  The  computations  were  done  on  a  CRAY 
Y-MP. 


Table  1.  Computation  Statistics  for  Three  Smooth  Surfaces. 

Smooth 

Maximum 

Average 

Number 

Computation 

Surface 

Deviation 

Deviation 

of 

Time 

(meters) 

(meters) 

hills 

(seconds) 

Ts(x,  y) 

5 

2.01 

630 

48 

Ux,y) 

2 

0.88 

3799 

236 

Ti(x,y) 

1 

0.34 

17,960 

1079 

If  T(x,y)  is  to  be  precomputed  and  then  used  in  a  simulauon,  the  difference  in 
computation  times  in  the  creation  of  the  smooth  surfaces  will  probably  not  be  too  impor¬ 
tant.  However,  since  the  computatiaial  burden  of  using  T(x,  y)  is  directly  proportional  to 
the  number  of  hills  comprising  the  surface,  "real  time"  applications  may  be  quite  sensitive 
to  this  factor. 

It  should  be  pointed  out  that  although  T(x,y)  was  constructed  using  the  VRT  hill 
format  (equation  2.4),  the  strict  seff-similarity  of  the  "ideal"  VRT  surface  was  not  used. 
While  we  are  free  to  determine  for  ourselves  the  distribution  of  terrain  features  when 
building  generic  terrain,  when  we  constmct  a  T{x,y)  to  match  a  given  set  of  DTAs,  the 
distribution  of  terrain  features  must  conform  to  the  peculiarities  of  that  particular  piece  of 
terrain.  Thus,  the  best  we  can  hope  for  is  approximate  self-similarity.  Even  this  may  not 
exist  if,  for  example,  the  set  of  DTAs  under  consideration  includes  two  or  more  different 
types  of  terrain. 

However,  when  it  comes  to  adding  micro-terrain  to  T{x,  y),  we  are,  in  effect,  creat¬ 
ing  a  miniature  set  of  generic  terrain.  Here,  the  full  power  of  the  ^'TlT  theory  can  be 
used.  In  general  terms,  the  empirical  distribution  of  the  VRT  parameters  comprising 
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Figure  3.  Sample  Contour  Map:  6=5 
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Figure  4.  Sample  Contour  Map: 
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T(x,y)  can  be  analyzed  and  used  to  randomly  generate  small  terrain  features,  whose 
scale  distribution  is  governed  by  equation  2.1. 


4.  Future  Work:  Optimization  of  the  Algorithm 

It  is  for  just  such  applications  that  optimization  of  the  algorithm  will  be  important. 
Here,  optimization  means  creating  T(x,  y)  (for  a  given  S)  using  the  fewest  hills  possible. 
Unfortunately,  for  a  given  set  of  DTAs,  we  will  not  know  what  that  minimum  number  of 
hills  is.  However,  there  are  some  techniques  we  can  apply  in  the  surface  construction 
process  that  will  tend  to  reduce  the  number  of  hills.  One  method  is  to  increase  the  num¬ 
ber  of  iterations  in  the  process  of  selecting  the  parameters  for  each  hill.  Early  experi¬ 
ments  indicate  that  this  does,  in  fact,  yield  a  reduction  in  the  number  of  hills,  although  it 
comes  with  an  increase  in  calculation  time. 

Another  idea  is  to  increase  the  variability  of  the  hill  parameters  in  the  iteration 
loop.  This  will  allow  us  to  explore  more  varied  hill  shapes  and,  presumably,  find  a  better 
fit  However,  the  rate  of  convergence  will  be  reduced,  making  it  necessary  to  increase  the 
number  of  iterations  per  hill  to  find  that  better  fit.  We  will  have  to  find  the  proper  balance 
between  increased  parameter  variability  and  number  of  iterations. 

It  is  also  possible  to  make  the  iterative  loop  process  adaptive.  By  this  we  mean  that 
the  rate  of  improvement  of  fit  during  the  iteration  loq?  will  be  monitored,  and  the  itera¬ 
tion  automatically  ter:ninated  when  the  rate  erf  improvement  falls  below  a  specified 
threshhold  value. 

Finally,  instead  of  treating  the  hill  fitting  process  locally  (i.c.,  by  fitting  the  hill  at 
the  site  of  single  highest  residual  DTA),  we  could  look  at  the  process  on  a  "regional" 
basis,  here,  we  would  do  a  global  statistical  analysis  of  the  residual  DTAs,  find  the 
largest  "ccrfierent"  topographical  feature,  and  try  to  find  a  hill  to  fit  the  whole  thing.  If 
this  could  be  done,  it  would  clearly  speed  up  the  entire  process.  The  open  question  is 
how  to  define  the  largest  "coherent"  topographical  feature. 

Also,  additional  woik  is  necessary  to  fully  develop  the  capabilit>  to  superimpose 
micro-terrain  on  T(x,  y). 
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RAY  NORVELLE 
GEORGE  LUKES 
7701  TELEGRAPHIC  ROAD 
ALEXANDRIA  VA  22310-3864 

1  INSTITUTE  FOR  DEFENSE  ANALYSIS 
ATTN:  WARREN  K.  OLSON 
ROOM  242W 

1801  N  BEAUREGARD  ST 
ALEXANDRIA  VA  22311-1772 


ABERDEEN  PROVING  GROUND 
58  DIR,  USARL 

ATTN:  AMSRUWT-W/CHARLES  H  MURPHY 
AMSRL-WT-WE/ 

JUDITH  TEMPERLEY 
MARY  ANNE  FIELDS 
CAROLYN  PATTERSON 
JOSEPH  WALD  (50  CP) 
AMSRL-a-SB/ 

VIRGINIA  KASTE 
THOMAS  KENDALL 
CHARLES  HANSEN 
TERRY  PURNELL 

2  DIR,  USAMSAA 

ATTN:  AMXSY-ED/IACK  MEREDITH 

AMXSY-CD/DWAYNE  NUZMAN 
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USER  EVALUATION  SHEET/CHANCE  OF  ADDRESS 


This  Lalwratory  undertakes  a  rx)nlinuing  effort  to  improve  the  quality  of  the  reports  it  publishes.  Your 
comments/answet^  to  the  items/questions  below  will  aid  us  in  our  efforts. 

1 .  ARL  Report  Number  ARL-TP-605 _ Date  of  Report  Ormhpr  1^94 _ 

2.  Date  Report  Received _ _ _ 

3.  Does  this  report  satisfy  a  need?  (Comment  on  purpose,  related  project,  or  other  area  of  interest  for 

which  Ute  report  will  be  used.) _ 


4.  Specifically,  how  is  the  report  being  used?  (Information  source,  design  data,  procedure,  source  of 
ideas,  etc.) _ 


5.  Has  the  information  in  this  report  led  to  any  quantitative  savings  as  far  as  man-hours  or  dollars  saved, 
operating  costs  avoided,  or  efficiencies  achieved,  etc?  If  so.  please  elaborate. _ 


6.  General  Comments.  What  do  you  think  should  be  changed  to  improve  future  reports?  (Indicate 
changes  to  organization,  technical  content,  format,  etc.) _ _ 


Organizatiem 


CURRENT  Name 

ADDRESS _ 

Street  or  P.O.  Box  No. 


City,  State,  Zip  Code 

7.  If  indicating  a  CJiange  of  Address  or  Address  Correction,  please  provide  the  Current  or  Correct  address 
above  and  the  Old  or  Incorrect  address  below. 


Organizatton 


OLD  Name 

ADDRESS _ 

Street  or  P.O.  Box  No 


City,  State,  Zip  Code 


(Remove  this  sheet,  fold  as  indicated,  tape  closed,  and  mail.) 
(DO  NOT  STAPLE) 


Department  of  the  Army 


OFFICIAL  BUSINESS 


BUSINESS  REPLY  MAIL 

FIRST  CLASS  PERMIT  NO  0001,  APG,  MD 


PoBtsgB  wiH  bB  paid  by  addmaaa 

Director 

U.S.  Army  Research  Laboratory 

ATTN:  AMSRL-OP-AP-L 

Aberdeen  Proving  Ground,  MD  21005-5066 


